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^ ■ Abstract 
CN ' 

If^ \ The spectral energy distribution (SED) of gamma-ray loud BL Lac objects typ- 

' icahy has a double-humped appearance usuahy interpreted in terms of synchrotron 

seh-Compton models. In proton blazar models, the SED is instead explained in 
Q I terms of acceleration of protons and subsequent cascading. We discuss a variation 

' of the Synchrotron Proton Blazar model, first proposed by Miicke Sz Protheroe 

• (1999), in which the low energy part of the SED is mainly synchrotron radiation 

by electrons co-accelerated with protons which produce the high energy part of the 



O ' SED mainly as proton synchrotron radiation. 

-4— > ' As an approximation, we assume non-relativistic shock acceleration which could 

, apply if the bulk of the plasma in the jet frame were non-relativistic. Our results may 

therefore change if a relativistic equation of state were used. We consider the case 
where the maximum energy of the accelerated protons is above the threshold for pion 
photoproduction interactions on the synchrotron photons of the low energy part of 
I the SED. Using a Monte Carlo/numerical technique to simulate the interactions and 

subsequent cascading of the accelerated protons, we are able to fit the high-energy 
gamma-ray portion of the observed SED of Markarian 501 during the April 1997 
flare. We find that the emerging cascade spectra initiated by gamma-rays from vr*^ 
decay and by from /j^ decay turn out to be relatively featureless. Synchrotron 
radiation produced by /i^ from vr^ decay, and even more importantly by protons, 
and subsequent synchrotron-pair cascading, is able to reproduce well the high energy 
part of the SED. For this fit we find that synchrotron radiation by protons dominates 
the TeV emission, pion photoproduction being less important with the consequence 
that we predict a lower neutrino flux than in other proton blazar models. 
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1 Introduction 



During its giant outburst in April 1997, the nearby BL Lac object Mkn 501 (at redshift 
z=0.034) emitted photons up to 24 TeV and 0.5 MeV in the 7-ray and X-ray bands, 
respectively, and has proved to be the most extreme TeV-blazar observed so far (e. g. 
Catanese et al 1997, Pian et al 1998, Protheroc et al 1998, Quinn et al 1999, Aharonian 
et al 1999). This energy is the highest so-far observed for any BL Lac object, and the 
flux is approximately 2 orders of magnitude higher than the synchrotron peak at its 
quiescent level. BeppoSAX and OSSE observations (Maraschi 1999) suggest that the X- 
ray spectrum is curved at all epochs, and the spectrum during fiaring has been fitted by a 
multiply-broken power-law (Bcdnarck & Protheroe 1999). COMPTEL has not seen any 
significant signal from Mkn 501 at any time (Collmar 1999), while a Sa upper limit of 
F{> 100 MeV) < 3.6 x lO^^cm-^ s"^ has been derived for the April 1997 EGRET viewing 
period (Catanese ct al 1997). 

A flux increase at TeV-energies was also observed with the Whipple, HEGRA and 
CAT telescopes (Catanese et al 1997), with the most intense flare peaking on April 16 at 
a level ~ 100 times higher than during its quiescent flux. The non-detection of Mkn 501 
by EGRET indicates that most of the power output of the high energy componcnt is in the 
GeV-TeV range. The TeV-observations revealed a power-law spectrum with photon index 
~ 2 up to ~ 10 TeV and a gradual steepening up to 24 TeV. The extragalactic diffuse 
infrared background Icads to significant extinction of 7-rays through 77-pair production 
above 10 TeV. The extinction-corrected TcV-spectrum (e. g. Bednarek & Protheroc 1999), 
shows the spcctral energy distribution (SED) peaking at ~ 2 TeV. Optical observations 
did not show any significant variations (Buckley & McEnery 1997), indicating that the 
change in the low energy part of the SED was mainly confined to the X-ray band above 
0.1 keV. 

Various models have been proposed to explain the observed 7-ray cmission from TeV- 
blazars, all of which are identified as high-freqiicncy pcakcd BL Lac objects. Leptonic 
models, in which electrons inverse-Compton scattcr a population of low energy photons 
to high energies, currently dominate the thinking of the scientific community. Because 
of the low luminosity of accretion disks in BL Lacs, the main target photons for the 
relativistic electrons would be the synchrotron photons produced by the same relativistic 
electron population, as in the synchrotron sclf-Compton (SSC) model. An alternative 
scenario for the production of the observed 7-ray fiux has been proposed involving pion 
photoproduetion by energetic protons with subsequent synchrotron-pair cascades initiated 
by decay produets (photons and e^) of the mesons (e. g. Mannheim et al 1991, Mannheim 
1993). Thcsc proton-initiated cascadc (PIC) models could, in principle, be distinguished 
by the observation of high energy ncutrinos produced as a result of photoproduetion. 

In this paper, we consider the April 1997 flare of Mkn 501 in the light of a modifled 
Synchrotron Proton Blazar (SPB) model. We assume that electrons (e~) and protons (p) 
are accelerated by Ist order Fermi acceleration at the same shock. The relativistic e~ 
radiate synchrotron photons which serve as the target radiation field for proton-photon 
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interactions, and for the subsequent pair-synchrotron cascade which develops as a result 
of photon-photon pair production. This cascade redistributes the photon power to lower 
energies where the photons escape from the emission region, or "blob," which moves 
relativisticaUy in a direction closely ahgned with our hne-of-sight. 

Until recently, this model was not able to reproduce the general features of the double- 
humped blazar spectral energy distribution (SED), but produced a rather featureless spec- 
trum (see e.g. Mannheim 1993), nor could it explain correlated X-ray/TeV-variability. 
Here, we present a comprehensive dcscription of our Monte-Carlo simulations of a sta- 
tionary SPB model, including all relevant emission processes, and show that this model is 
indeed capable of reproducing a double-humped SED as observed. Here, the origin of the 
TeV-photons are proton synchrotron radiation, as first proposed by Miicke & Protheroe 
(1999); a similar model has also bccn proposed by Aharonian (2000), and Rachen (1999) 
presented speculations about /i^- and proton-synchrotron radiation Icading to narrow 
cascade spectra during fiares, which might explain correlated X-ray/TeV-variability Jet 
energetics and limits from particle shock acceleration, however, put severe constraints on 
this scenario. The goal of this paper is to discuss the physical processes included in our 
SPB model Monte-Carlo code, and givc the results of applying this code, as an example, 
to reproduce the SED of the giant fiarc from Mkn 501 which occurred in April 1997. A 
comprehensive study of the whole parameter-space (magnetic field, Doppler factor, etc.) 
for this model will be the subject of a subsequent paper. 

In Section 2, we discuss constraints on the maximum particle energies imposed by the 
co-acceleration scenario, and by the pion production threshold. Section 3 is dcvoted to 
the emission processes in the present model. Energy losses and particle production are 
treated in Sect. 3.1, while the cascade calculations, including a brief description of our 
code, are outhned in Sect. 3.2. In Sect. 4 we apply our model to the April 1997 fiare 
of Mkn 501. The multifrequency photon spectrum is shown in Sect. 4.1, while in Sect. 
4.2 the predicted neutrino spectrum is discussed. We conclude with a discussion and 
summary in Section 5. 

2 The Co-acceleration Scenario 

In the present model, shock accelerated protons (p) interact in the synchrotron photon 
field generated by the electrons (e~) co-accelerated at the same shock. This scenario may 
put constraints on the maximum achievable particle energies. 

The usual process considered for accelerating charged particles in high energy as- 
trophysics is diffusive shock acceleration (see e.g. Bell 1978, Drury 1983, Blandford & 
Eichler 1987, Biermann & Strittmatter 1987, Jokipii 1987, Jones & EUison 1991), in 
which particles undergo coUisionless scattering, e.g. by Alfven waves, in the upstream 
and downstream plasma. Charged particles with gyroradii larger than the thickness of 
the shock front propagate with diffusion coefficients Ki and k,2 in the upstream and down- 
stream plasma, rcspectively, for propagation parallel to the shock normal. In the shock 
frame the plasma fiow velocity changes from ui = j3ic in the upstream region to U2 in the 
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downstream region. In this paper, for simplicity we restrict ourselves to non-relativistic 
shocks, and postpone discussion of relativistic shock acccleration to a later paper. 
The acceleration time scale for non-relativistic shocks is given by 

'-^ = ^;r^^i^^ + r,^^) (1) 

where Tc — > 4 is the compression ratio for the case of strong shocks in a non-relativistic 
monoatomic ideal gas. If the magnetic field is governed by an ordered component, the 
orientation of the shock normal to the main magnetic field direction becomes important. 
In general, the diffusion coefficient can be written as 



K/'i 



^ Ki^llC0S^9i + Ki^j_sm'^ 9i, i = 1,2 (2) 



where 9i is the angle between the magnetic field and the axis connccting the upstrcam 
(i = 1) and downstream [i = 2) regions. In the diffusion limit, kinetic theory relates the 
parallel and perpendicular diffusion coefficients through 

i^W ^ [i + v'^]i^± ^ ^h^c (3) 

where 77 = \\/rg with A|| bcing the mcan frcc path parallel to the magnetic field, Vg = 
(3'ymc^/eB is the particle's gyroradius, m and 7 = (1 — are the particle's mass 

and Lorentz factor, respectively, and B is the magnetic field strength in the upstream 
region. The mean free path is, in general, a function of the particle energy through its 
gyroradius, and is dependent on the spectrum of the magnetic turbulence. In the small 
angle scattering approximation (i. c., if Alfven waves dominate the particle defiection with 
wavelength equal to the particle gyroradius; see Drury 1983) we have 

A = g (A\ 

8nl{k)k' ^ ' 

This spectrum 7(/c) is usually expressed as a power law of the wave number k in the 
turbulent magnetic field: 

I{k) oc k-^. (5) 

5 = 5/3 corresponds to Kolmogorov turbulence which may be common in astrophysical 
environments (Biermann & Strittmattcr 1987), while 5 = 1 corresponds to Bohm diffusion, 
and is often considered for simplicity. For strong magnetic fields, Kraichnan turbulence 
5 = 3/2 (Kraichnan 1965) may be present. In the foUowing, we consider 5 as a free 
parameter. The mean free path may then be expressed as (see Biermann & Strittmatter 
1987) 



<5-l 

.ax I 

6(5 - 1) V r~ 



A,| = -/^ r:^! for 5^1 (6) 
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where b is the ratio of the turbulent to ambient magnetic energy density, rg,inax is the 
gyroradius of the most energetic protons, and has the same order of magnitude as the 
system size, and r^^min corresponds to the smallest turbulence scale. The mean free path, 
and consequently the acceleration time at maximum energy is only shghtly dependent 
on the turbulence spectrum in the case of protons, whereas for electrons the acceleration 
time at maximum energy shows a strong dependence on the magnetic turbulence spectrum 
adopted. 

We expect 6 < 1, since otherwise the energy density in particles would not bc able to 
be confined by the ambient field (Biermann & Strittmatter 1987). With these relations, 
the acceleration time scale may be re-written as 



where 



rc-1 



cos^ Oi 



sin^ r, COS-' 0^ + siir' e j{l + z/^) 



(l + r?2) (cos2^i+r2sin2^i)3/2 



(7) 



(8) 



The diffusion approximation used here limits the maximum mean free path io r} < j3/ j3i 
(Jokipii 1987). If the particle spectra are cut off due to synchrotron losses, balancing the 
acceleration time scale with the loss time scale determines the maximum Lorentz factors 
of protons, 

7p,max = 2.1 X lO^i/?! [l3BF{eur]p,u.^.)]~^'^ (9) 

and electrons 

7e,max = 1-2 X 10^1 [l3 B F{e^,7]„)]-^''' (10) 

where B is in Gauss. The ratio of the maximum proton Lorentz factor to the maximum 
electron Lorentz factor is then 



< 



rrip [" F(6'i,?7e,max) 
rrie |.-F(6'l,77p,max) 



1/2 



where the equality corresponds to the maximum proton energy being determined by syn- 
chrotron losses, and the inequality to the maximum proton energy being determined 
instcad by adiabatic losses. For parallcl shocks this relation is consistent with the results 
found by Biermann & Strittmatter (1987). 

The corrcsponding acceleration time scales at the maximum particle energies, if de- 
termined by synchrotron losses, are 



t 



acc,p,max 



2.2xlO^S-3/2/3-i[/?F(e,,ry^^ 



,1/2 



(12) 
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for protons, and 

tacc,e,max = 6.68-'^^ f3^' , r/e,max] S (13) 

for electrons. In the present paper, we shall adopt jSi = 0.5, (3 = 1 and Vc = 4. The 
geometry dependent term [-^(6*1, ?7p^niax)]^^^ is plotted in Fig. la as the dashed curves for 
different ?7p,max values and Vc = 4. It can be seen that highly obhque shocks allow proton 
acceleration on very short time scales. In the hmiting case of a perpendicular shock, drift- 
shock acceleration drives the energy gain, and the finite size of the shock front restricts 
the maximum particle energy. The maximum proton Lorentz factor 7p is reached for the 
maximum drift distance, the shock size, and is given by 

7p,max = 3.2 X lO-'PiRB (14) 

with R in cm and B in Gauss. 

At their maximum Lorentz factors, the acceleration process for protons is considerably 
slower than for electrons, and so the proton acceleration time must be consistent with the 
observed variabihty time, tvar, 

^var-0 ^ ^acc,p,max (-^^) 

where D the Doppler factor. This can be converted to a constraint on the geometry 
dependent term using Eq. |12|, and is plotted as a function of 6i in Fig. la for a typical 
set of TeV-blazar parameters (i? f» 20 G, D ~ 10, tvar = 12 hours). The region below the 
solid line is allowed by the variabihty constraint (Eq. |T5|), and gives for each 77p,max value 
a minimum shock angle between 69° and 88°, depending on ?7p_max- Thus, proton shock 
acceleration on hour time scale in hadronic models can only take place in oblique shocks, 
and the maximum drift distance may restrict the maximum energy gain rather than the 
gyroradius. 

Using the same shock geometry and magnetic turbulence spectra for both protons and 
electrons, we find that the ratio -^(6*1, ?7e,max)/-P'(6'i, ^p,max) does not vary by more than a 
factor of 2 within the allowed shock angle range (see Fig. Ib), allowing us to adopt an 
average value for this ratio for a given parameter combination. Inserting this ratio into 
Eq. |Tl| then restricts the ratio of the allowed maximum particle energies to the range below 
the solid lines shown in Fig. 2. Points exactly on this line represent models where the 
particle spectra are limited by synchrotron losses and the acceleration time scale is exactly 
the variabihty time scale; points below this line apply if adiabatic losses are dominant for 
protons or the proton acceleration time scale is shorter than the variabihty time scale in 
the jet frame. We note that for a given maximum proton energy, the highest maximum 
electron energy occurs with Bohm diffusion. 

In hadronic blazar models pion photoproduction is essential for neutrino production. 
The threshold for this process is given by emax7p,max = 0.0745 GeV where emax is the 
maximum photon energy of the target field, which in the Synchrotron Proton Blazar 
models is in turn produced by the co-accelerated e~. In the (5-function approximation 
for the synchrotron emission, emax = {'i/^)l1m&AB / Bcv)TneC^ with Bcr = 4.414 x 10^^ G. 
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Inserting einax into the threshold condition, we find 



,„_>1.72xl0-(^)"7.-„ (16) 



and this is shown in Fig. 2 as the dashed lines for various magnetic field strengths. To- 
gether with Eq. 11 the allowed range of maximum particle energies is then restricted to 
the area below the sohd hnes and abovc the dashed hncs in Fig. 2, as shown for cxam- 
ple as the shaded area for B ^ 20 G and Kolmogorov turbulence. Inspecting hadronic 
models as presented in the hterature (e. g. Mannheim 1993, Mannheim et al 1996, Rachen 
1999), we find that most models which are able to fit the observations he above the Bohm 
diffusion hne (see Fig. 2), indicating that the turbulence spectrum required in common 
hadronic blazar jet models is likely to be of Kolmogorov/Kraichnan type. Rachen (1999) 
speculates that the transition bctween Kolmogorov and Kraichnan type turbulence could 
be responsible for the difference between low- and high frequency peaked BL Lacs. 



3 Emission Processes 

In the present model, the co-accelerated e~ are assumed to produce most of the observed 
low energy part of the SED by synchrotron radiation, and this is assumed to be the target 
radiation field for p7 interactions and subsequent cascading. The observed hardening of 
the spectrum with rising flux, has recently been convincingly rcproduced by a shock model 
with escape and synchrotron losses (Kirk et al 1998). The spectral slopes in this model 
are controUed by synchrotron cooling, and thus naturally explain the temporal behaviour 
of the spectral index as observed in the X-ray band. The flaring behaviour is explained 
by the shock front running into plasma whose density is locally enhanced, and which 
thus increascs the number of particlcs injected into the acceleration proccss. The increase 
in the plasma density is accompanied by an increase in the magnetic field, leading to a 
higher acceleration rate and a shift of the maximum particle energies to higher energies. 
This picture departs from the Standard explanation in which tlic flat synchrotron spectra 
appear as a result of a superposition of several local self-absorbed synchrotron spectra 
with changing self-absorption frequency, adopted in previous PIC models (e.g. Mannheim 
1993). 

For simplicity, and because we do not wish to include additional parameters, we use 
the same magnetic field for synchrotron radiation as for acceleration. For normal shocks 
this approximation is justified as the magnetic fields either side of the shock are similar. 
This approximation might evcn be justified for oblique shocks as a lower magnetic field 
in the upstream region, compared to the downstream region, implies a higher diffusion 
coefficient and time spent upstream, increasing the synchrotron losses there and partially 
compensating for having a lower field upstream. Also, at oblique shocks, refiection at 
the shock front itself is thought to be more important than diffusion in the downstream 
region (Kirk & Heavens 1989), so that accelerating particles spend most of their time 
upstream. In addition, we assume that pitch-angle scattering maintains quasi-isotropic 



7 



particle distributions, and all radiating particles are confined to the homogeneous emission 
region. 



3.1 Energy Losses 

There are several energy loss/interaction processes which are important for protons, 
electrons and photons in a dense radiation field produced by relativistic electrons co- 
accelerated along with the protons: protons interact with photons, resulting in pion pro- 
duction and (Bethe-Heitler) pair production; electrons, muons, protons and charged pions 
emit synchrotron radiation; photons interact with photons by pair production. We shall 
show below that, for the present model, Inverse Compton emission by the electrons can 
be neglected. 

For simplicity we represent the observed synchrotron spectrum of Mrk 501 during 
haring, the target photon field for the p7-collisions and photon-photon pair production, 
as a broken power-law: 

r for < e < 1.6 kcV , . 

H^) « I g-1.8 fo^ 1 g keV < e < 42 keV 

For determining the photon density of the target field, the dimcnsion of the emission 
region, assumed to be spherical, must be known. This can be estimated by setting the 
photon Crossing time equal to the variability time scale (in the jet frame - in the remainder 
of this section all quantities are in the jet frame unless noted otherwise), making the 
implicit assumption that the light crossing time scale determines the fiux variations. The 
observed variability time scale can, in general, depend on: (i) the injection time scale 
for the energetic particles tacc(-E); (ii) the time needed for converting their energy into 
radiation, i. e. their energy loss time scale tioss{E); (iii) the effcctivc light crossing time 
tcross(-E') ~ (2-Rbiob/c) X Pescl-^') where i?biob is the geometrical blob radius and Pesc{E) is 
the energy dependent probability for photons escaping from the blob taking account of 
77-pair production and diffusion during cascading. Hence, 

-^^var ~ I^3'X(^acc) ^loss) ^cross) ■ (-'-S) 

For leptonic models, the time scales for energy losses and acceleration are typically 
significantly shorter than the crossing time, i. e. tiossj^acc ^ ^crossj and thus the radius of 
the emission region can be derived from the observed variability time scale. In addition, 
the emission region is assumed to be optically thin at 1 TeV and X-ray energies, implying 
that i?TeV ~ Rx, where R^eV, Rx are the dimensions of the emitting region at 1 TeV and 
X-ray energies, respectively. This differs in two points from the hadronic blazar jet models: 
Firstly, the optical dcpth of the emission region is strongly energy dependent, leading to 
an effective, i. e. observed, thickness of the emission region Res{E) ^ i?biob-Pesc (-£'), which 
also depends on the energy. Diffusion can be neglected during cascading since the cascade 
processes are of leptonic origin, and are, in general, more rapid than diffusion. Thus, in 
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SPB models the crossing time scale in the optically thick TeV-band is related to the 
Crossing time scale at X-ray energies (optically thin), through 

^cross,TeV ~ ^cross,x[l " exp(— r^^^Tcv)]/T77,TeV (19) 

(avcraging over a homogeneous emission volume; see Rachen 1999). Hcncc, tcross,x is the 
relevant time scale for estimating the radius -Rbiob ~ Rx of the emission region. Sccondly, 
the acceleration and/or energy loss time scales can be of the same order of magnitude as 
the crossing time scale. Acceleration, however, is always faster than the energy losses of 
the accelerating particles up to their maximum energy. Note that because of the leptonic 
nature of the cascadc proccsscs, the energy loss time scale will in general be determined 
by the (slowcr) hadronic proccsses. As a further consequence, fiux variations are not 
significantly washed out by the cascading mechanism, but closely foUow the crossing or 
hadronic loss time scales. 

If ^loss ^ ^cross the crossing time scale determines the fiux variations, and in this case we 
can estimate the radius of the emission region through -Rbiob ~ 0.5cDtvar,x with tyar.TcV ~ 
^var,x/T77,TcV- OUT prcscut modcl, efficient proton synchrotron emission determines the 
TeV-bump in the blazar SED, and so at the maximum proton energy tad > ^p,syn- The 
adiabatic loss time scale due to expansion is (Longair 1994) 

= = |i?biob/i?biob| - 2|5/^| « R/u, (20) 

assuming magnetic fiux conservation B oc R^^^,, and so t^d is related to the size of the 
emission region -Rbiob- Sincc ti^ss ~ *p,syn for the highest proton cncrgics in our model 
responsible for the TeV-cmission, tioss < ^cross) ^-^d so the size of the emission region can 
be determined by the variability time scale. 

The shortest doubling time measured by the Whipple Telescope in the 1997 data of 
Mkn 501 was approximately 2 hours (Quinn et al 1999), while HEGRA reported a lower 
doubling time of 15 hours (Aharonian et al 1999) or 12 hours (Krawczynski 1999). As a 
working hypothesis we adopt here a variability time scale of 12 hours, but will discuss also 
effects of smaller variability time scales. For tyar ~ 12 hours we find i^biob ~ 8 x lO^^cm 
for D ^ 10. 

Eq. 12 and 13 imply that the acceleration time scale for electrons is much smaller 
than for protons at their maximum energies, tacc,e,max <S ^occ,p,max- If a 7-ray outburst 
thus corresponds to particle acceleration of a singlc p/e~ population, then an increase of 
the low energy synchrotron fiux will occur before the TcV-fiarc. The amount of the lag 
depends on the proton acceleration time, and so on the shock parameters, in particular 
the diffusion coefficients. E. g. for the present model and the Kolmogorov spectrum 
tacc,p,max ~ 10^ S (scc Fig. 3), which for D=10 suggests a time lag of f» 5 hours. For Bohm 
diffusion, the time lag would be ~ 14 hours. These time lags are consistent with current 
observational constraints. 

The relevant radiation and loss time scales for photomeson production, Bethe-Heitler 
pair production, p synchrotron radiation, and adiabatic losses due to jet expansion, are 
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shown in Fig. 3 together with the acceleration time scale. Synchrotron losses, which 
turn out to bc at Icast as important as losses duc to photopion production in our model, 
limit the injcctcd p spcctrum to a Lorentz factor of 7p 3 x 10^° for the assumed 
model parameters. We adopt a Kolmogorov spectrum of turbulence for the magnetic 
field structure {S — 5/3), and so for any rjp > 1 value, variability arguments constrain 
the shock angle to 9i > 75° (see Fig. 1). The maximum proton energy could then be 
achieved, e.g., with rjp — 10, 6i = 85° and ui = 0.5c. This is in agrccment with the limit 
imposed on quasi-perpendicular shocks due to their finite shock size. Note that due to the 
non-zero shock angle, the acceleration time scale shown in Fig. 3 does not follow a strict 
power-law, but is curved. This is due to the non-linear dependence of F{9i,r]p) on the 
particle's gyroradius (see Eqs. 6-8). A 7"^ proton spectrum, typical of shock accelerated 
particles, is used for 2 < 7^ < 7p,max, where 7p,max = 3 x 10^° is obtained from requiring 

Rachen & Meszaros (1998) noted the importance of synchrotron losses of /j,^ (and tt^) 
prior to their decay in AGN jets and GRBs. The critical Lorentz factors 7^ ?a 2 x 10^ and 
7,^ 4 X 10^°, above which synchrotron losses in the assumed magnetic field dominate 
over decay, lie below the maximum Lorentz factor for /j,"^ and above the maximum Lorentz 
factor for tt^. Thus, while vr^-synchrotron losses can safely be ignored, /^''^-synchrotron 
losses should be included. 

For the parameters employed in this work we receive a target photon energy density 
of i^target ~ 0.06 TeV/cm^, and a magnetic field energy density of ub ~ 11-7 TeV/cm^ 
using B ^ 20 G. With ub ?^target significant Inverse Compton radiation from the 
co-accelerated er is not expected. 

We also neglect secondary production interactions of relativistic protons with the 
ambient thermal plasma. In the present model the dominant loss process turns out to 
be proton synchrotron radiation. Thus, our assumption is justified if the thermal proton 
density does not exceed k, lo'^- -^'' cm~^. In Section 4 we estimate the number density of 
cold protons to be less than this for reasonable values of the jet width. 

3.2 Simulation of particle production and cascade development 

For the first time in the context of the SPB-model, we use the Monte-Carlo technique 
to simulate particle production and cascade development, and this allows us to use exact 
cross sections. For photomeson production we use the Monte-Carlo code SOPHIA (Miicke 
et al 2000), and Bethe-Heitler pair production is simulated using the code of Protheroe 
& Johnson (1996). We calculate the yields for both processes separately, and the results 
are then combined according to their relative interaction rates. 

The mean pion production interaction rate for an isotropic photon field is 

'^AEp) = / de—r / dsp{sp-mp)a^{sp) , (21) 

where Sp = + 2Epe{l — Pp cos9p) is the center-of- momentum (CM) energy squared, 9p 
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the angle between the proton and the photon, f3pC the proton's velocity, ^/Sp^th ~ l.OSGeV 
the threshold CM energy, eth = (sp,th - ml)/2{Ep + pp), Sp,max = ml + 2Epe{l + (3p) and 
cr^ the pion production cross section. The Bethe-Heitler pair production interaction rate, 
Rbh, is calculated using the formulae given in Chodorowski (1992). 

In highly magnetized environments, proton-photon interactions compete with syn- 
chrotron radiation by the protons. To take proton synchrotron losses into account in 
our code, we sample a pj interaction length from an exponential distribution with its 
corresponding mean Xp^ = c/tt^^bh as given below. We make the approximation that the 
proton energy on interacting is given by 



E ^ 



■'^ p ,syn\^p j 



(22) 



where E^^ is the initial proton energy, Xp^ is the sampled interaction length for pion pro- 
duction or Bethe-Heitler pair production, and Xp,syn = c/rp^syn is the proton synchrotron 
loss distance given by 

^M,s,n = -[j^) -j^. (23) 

where M = rup, is the Thomson cross section, and ub is the magnetic energy density. 

The relatively long mean life time of highly energetic charged pions, muons and kaons 
might be of the same order of magnitude as their synchrotron loss time scale in highly 
magnetized environments like AGN jets and GRBs (Rachen & Meszaros 1998). We sim- 
ulate their synchrotron energy losses by sampling the decay length Xdec from an exponen- 
tial distribution with corresponding mean decay length Sdec = C'-yr^ec- For jji^ we have 
Tdec = 2.20 X 10-*^ s, while r^ec = 2.60 x 10"^ s for 7r± and r^ec = 1-24 x 10"^ s for K^. 



The particle's energy on decaying is then 



2x 



dec 



^K,n,ii,syn{Ej^ Ij) 



-1/2 

(24) 



where E^'^^^^^ is the initial K,7r, or /i energy, XK,n,fi,syn = c/ r k^-k, ^,syn is the synchrotron loss 
distance with rK,TT,^i,syn given by Eq. ^ for M = rriK, m^, or m^. 

We next outline the simulations of the cascade development in more detail. Energetic 
photons will pair produce on the target photon field, and if the magnetic energy density 
exceeds the injected target field density (mb > t^targct) they will initiate a pair-synchrotron 
cascade. Energetic photons arise directly from 7r°-decay, or indirectly as synchrotron 
photons from protons, charged mesons and electrons resulting from vr^ — /i^ — decay 
and Bethe-Heitler pair production. The optical depth t^^{E^) for 7-ray photons with 
energy E-^ for pair production inside the blob is given by 

^77(^7) = / c?e— ^ / ds^s^(T^^[s.,) (25) 

^-^/y ^ ^min ^ -^7, min 
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where n{e) is the differential photon number density and c-y-y(s) is the total cross section 
for photon-photon pair production (Jauch & Rohrlich 1955) for a centre of momentum 
frame energy squared given by 



(26) 



where 6y is the angle between directions of the energetic photon and soft photon, and 
S7,min = (2meC^)^, S7,max = 4ei?^ and emin = ijneC^Y / E^. For simulating photon-photon 
pair production we approximatc Ee+ = E^- = E^/2. The radiating are assumed to 
be continuously isotropized in the blob frame by deflection in the uniform magnetic field, 
resulting in the synchrotron radiation being isotropic in the blob frame. The spectrum 
of synchrotron photons, averaged over pitch- angle because of the isotropy of the particle 
distribution, is calculated using functions givcn by Protheroe (1990). 

In general, the cascade can be initiated by photons from vr^'-decay ("vr^ cascade"), 
electrons from the tt^ — > ^ decay ("tt''^ cascade"), from the proton-photon 
Bethe-Heitler pair production ( "Bethe-Heitler cascade") and p and //.-synchrotron pho- 
tons ( "p-synchrotron cascade" and "//^-synchrotron cascade"). Here, we assume the 
cascades develop linearly, and this requires the photon ficld produced by the cascade, 
to be negligible as a target field in comparison with the injected synchrotron radiation 
due to co-accelerated e~. This requirement can be expressed as r-y-y^cas 'C r-^^^^target with 
'r77,target)T77,cas bciug the pair productiou optical depths of photons on the target and cas- 
cade photon field, respectively. Our Monte-Carlo results show this condition to be met 
for the present input. To simplify the calculation, the electrons are completely coolcd in- 
stantly by synchrotron radiation before pair production by the synchrotron photons takes 
place. This approximation is equivalent to assuming that tgyn <^ tpair which is justified 
because of the very short synchrotron life time of electrons in the assumed magnetic field. 

A matrix method (e.g. Johnson et al 1996) is then used to foUow the pair-synchrotron 
cascade in the ambient synchrotron radiation field and magnetic field. The cascades are 
considered in the jet frame. Here, electron and photon fiuxes are represented by vectors 

and , which givc the total number of electrons in the energy bin at energy Ej and 
number of photons at energy Ei, respectively, in the /cth cascade gencration. We use a 
logarithmic stepsize of 0.1 ranging from log(£'e/l GeV) = —3 to 12 in electron energy, and 
from log(£'^/l GeV) = —13 to 12 in photon energy. Averaged over a homogeneous emis- 
sion region, the probability of gamma-ray interaction by photon-photon pair production 
at energy E^ is given by the vector P^^ i 



The transfer matrix T^^ gives the number of synchrotron photons of energy Ei 
produced by electrons of energy Ej^ Tj;j^^\ and the transfer matrix T^J^^^^ gives the number 
of of energy Ej produced through pair production of energetic photons of energy 
Ei on the target field. The vectors and matrices are calculated taking care of energy 
conservation. 



P. 



[1 - Pesc(^»)] = 1 



1 - exp(-r^^(-Ej)) 

Xy^(£'j) 



(27) 
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The photon and electron fluxes due to synchrotron radiation and photon-photon pair 
production are then calculated through matrix multiplication. We start by calculating 
the number of 7r°-decay photons which pair produce in the blob 

where the vector J? gives the number of 7r°-decay gamma-rays at an energy (the 
emerging photons, i. e. those which do not pair produce in the blob, are stored in an 
array). The electron spectrum due to photon-photon pair production is then given by 

^1 _ ^OjnCpair) 
j Z-^ i i j 

i 

These electrons radiate synchrotron photons, and the resulting photon yield is determined 

by 

7-1 _ \ ^ T^lrp(syn) 

j 

This is the Ist cascade generation photon spectrum, which in turn again suffcrs photon- 
photon pair production on the target photon field in the blob, etc. In our calculation, we 
iterate until F^^ < 0.01 at any energy Ei, or stop at 10 generations. 

We have injected 10^ protons at each proton energy equally spaced in logE'p at 0.1 
decade intervals from 10^ GeV to the maximum injection energy. These protons are as- 
sumed to be continuously isotropized, and interact with the synchrotron radiation field 
of the co-accelerated e~ through pion production and Bethe-Heitler pair production. The 
resulting neutrinos (see Sect. 4.1) escape without further interaction, but the and 
high energy 7-rays initiate cascades. The resulting particle spectra were weighted with 
Up^oEp'^dEp, appropriate to the assumed proton injection spectrum, and divided by the 
number of injected protons, and by An steradians. Figs. 4 and 5 show cxamples of cas- 
cade spectra in the observer's frame initiated by photons with different origins for model 
parameters (given above), which satisfactorily reproduce the fiare spectrum of Mkn 501. 

Our Monte-Carlo program only treats the first interactions of protons, and so we 
must take account of this when making the final fiux predictions. In the case of pion 
production, a fraction k ~ 0.25 of the initial proton energy gocs into particle production, 
and ~ 1/3 of the time the emerging nucleon is a neutron, and is assumed to escape from 
the blob without further interaction. Combining these two factors, the resulting spectra 
must be multiplied by 2.01 to take account of subsequent pion production interactions. 
Similarly, the multiplication factor which takes account of subsequent Bethe-Heitlcr pair 
production interactions is ~ 10. Particle spectra due to pion production and Bcthe- 
Heitler pair production are then weighted according to their mean interaction rates, i. e. 
-Rtt/ Rtot for pion production, and Rbh/ Rtot for Bethe-Heitler pair production where i?tot — 
Rt^ -\- it^BH + Rad is the total interaction rate, taking into account adiabatic losses of the 
protons due to jet cxpansion approximated by Eq. 18. 

Adiabatic jet expansion also affects the radius of the emission region, its photon dcnsity 
and magnetic field, and thereby the probability for p7-interactions and the subsequent 
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cascade development. With the expansion of the emission region, the probability for p7- 
interactions decreases because of the decrease of the photon energy density and magnetic 
field. It foUows that p7-interactions, and their resulting cascades, most hkely take place 
during the initial phase of the jet expansion, and so we neglect any effects on the emerging 
cascade spectra due to changes in magnetic field or dimension of the emission region. In 
particular, the fiare rise-time can be much shorter than the expansion and hght crossing 
time scales if one associates it with the onset of the p7-interactions and their resulting 
cascades. Short rise-times and longer flare decay time scales are typically observed in 
blazar light curves. Finally, we transform the cascade spectrum to the observer's frame. 
For calculating the distance, qo = 0.5 and Ho = 50 km s~^ Mpc~^ are used. 

Cascades initiated by photons from 7r° decay, "vr^ cascades", or by electrons from 
vr — > |U — i> e decay, "7?=*= cascades", (Fig. 4) produce rather featureless spectra (see also 
Mannheim et al 1993). However, cascades initiated by protons "p-synchrotron cascades", 
and by muons from vr — /i — > e decay, "/i^-synchrotron cascades" (Fig. 5b) produce a 
double-humped SED as observed for 7-ray blazars (see also Rachen 1999). The contribu- 
tion from Bethe-Heitler pair production turns out to be negligible (Fig. 5a). In the model 
for Mkn 501 presented here, the synchrotron flux from directly accelerated e~ significantly 
exceeds the reprocessed /x and p synchrotron flux. Thus, in the present Mkn 501 model 
direct proton and muon synchrotron radiation is mainly responsible for the high energy 
hump whereas the low energy hump is mainly synchrotron radiation by the directly ac- 
celerated e~ . In the present model photons up to 1 TeV in the jet frame are optically 
thin to 77-pair production. Thus, only a small fraction of the power emitted as proton 
synchrotron radiation is redistributed to lower energies via cascading. For models where 
the optical depth of TeV-photons is significantly higher, a correspondingly larger fraction 
of the TeV-bump may be redistributed to X-ray energies. In these models, the pairs 
produced by photons of the "high energy hump" may contribute considerably to the ob- 
served X-ray bump, whereas in the present model they are negligible in comparison to the 
directly accelerated electrons. We note that for the case of strong magnetic fields where 
/i='=-synchrotron cascades cause detectable humps, the p cascade spectrum dominates in 
general over the /^^ cascade spectrum. 

4 Application to the April 1997 flare of Mkn 501 

Adding the four components of the cascade spectra in Fig. 4-5 we obtain the SED shown 
in Fig. 6 where it is compared with the multifrequency observations of the 16 April 
1997 flare of Mkn 501. We use the parametrization of Bednarek & Protheroe (1999) to 
represent the BeppoSAX+OSSE data (thick straight line at low energies). The broken 
power-law simplification (Eq. 0) used as the target radiation field for p7-interactions and 
77 interactions is shown by the chain line. The 100 MeV upper limit from Catanese et al 
(1997) is nearly simultaneous to the TeV-fiare observations. The TeV-fiux, corrected for 
pair production on the cosmic background radiation field (Bednarek & Protheroe 1999) 
for two different background models, is shown as the thick curves. 
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The parameters used for modeling the April 1997 flare are: D — 12, B ki 20 G, radius 
of the emission region i?biob = 8 x 10^^ cm. For a target photon ficld for the p7-interactions, 
and the cascades as given in Eq. 16, wc find a photon energy density of this radiation field 
of «target = 60 GeV/cm~^. The accelerated protons are assumed to follow a power law 
oc between 2 < 7^ < 7p,max = 3 x lO"*^*^, and in order to fit the emerging cascade spectra 
to the data a proton number density of Up ~ 250 cm~^, corresponding to an energy density 
of accelerated protons of Up ~ 11.6 TeV/cm~^ is required. With a magnetic field energy 
density of Ub ~ 11.7 TeV/cm~^ our model satisfies «target ^ Up ^ Ub (all parameters 
are in the co-moving frame of the jet), confirming that a significant contribution from 
inverse-Compton scattering is not expected. 

To calculate the total jet luminosity Ljet, measured in the rest frame of the galaxy, 
we adapt the formulae of Bicknell (1994) and Bicknell & Dopita (1997) given for the 
synchrotron self-Compton model to apply for the case of our proton blazar model. We 
receive 



T-high 
-^obs 

DKp 



(r - 1) PB CpS'z 



(28) 



where r = (1— /3^)~^ ^ D/2\s a good approximation to the Lorentz factor of jets closely 
aligned to the line of sight, and the four terms inside the bracket give the relative con- 
tributions to the total jet power of cold protons, accelerated protons, magnetic field, and 
accelerated electrons, respectively. The contribution from cold protons is estimated on the 
basis of charge neutrality. 5'^'^" and S^^^^^ are the observed bolometric fiuxes of the low and 
high energy component, respectively, and Ce ~ 1, Cp ^"^^ the radiative efficiencies for elec- 
trons and protons. iJ^^^ — Ai^dFS^^^ with d the luminosity distance, ps — [-B^/(2//o)]/3 
and 



L 



high 

obs 



gives the jet-frame pressure of relativistic protons that would apply in the absense of 
energy loss mechanisms, and 

^ 3 / TUp Cp'S'ob" 1 1 A .„ s 

4UeCe5„^^^7elln(7e2/7el) l,MlP2hpi) ) ' 

We use this formula to estimate the total jet power to be ~ lO''^ erg/s, and we find 
the contributions to the total jet power of cold protons, magnetic field, and accelerated 
electrons, relative to that of accelerated protons, to be 0.8, 1, and 0.01, respectively, with 
a number density of cold protons ~ 10^ cm~^. Thus, for the shock acceleration mechanism 
we require more power per particle going into relativistic protons than into electrons. 

Fig. 7 shows the dependence of the total jet luminosity on the Doppler factor D with 
fixed parameters tyar = 12 hours, (3i = 0.5, and with B, R, Up and «target being of Mkn 501 
during fiaring. Clearly visible is the fact that in hadronic models the proton kinetic energy 
and the poynting fiux dominate the total jet luminosity, while the electron kinetic energy 
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is only of minor importance. At high Doppler factors the emission region becomes so 
large that one nccds only rclativcly small magnetic fields and proton densities to fit the 
observations. In addition, adiabatic losses become small resulting in a decrease of the 
required kinetic proton luminosity. For example, 5 5 G and Up ~ 10~^ cm~^ are 
sufficient to fit the Mkn 501 fiare for D — 50, while for D — 8 magnetic fields of over 
30 G and proton densities of Up f» 10^ cm~^ are needed. The total jet luminosity exhibits 
a minimum of ~ lO^^erg/s at around D Ki 12, which corresponds to the model we have 
chosen to present here. 

Smaller variability time scales in the shock acceleration model can only be explained 
by quasi-perpendicular shocks with a large mean free path (see Fig. la), which in turn 
are only consistent with the diffusion approximation for low shock speeds. Due to the size 
limit of the shock, the maximum energy gain is then limited to significantly lower values, 
and one may not reach TeV-photon energies in the proton synchrotron model, unless the 
magnetic field is increased accordingly, which in turn leads to larger jet luminosities. For 
example, for D — 10 and a variability time scale of 12 hours, at least B=20 G is needed 
to comply with all acceleration constraints. For t^ar = 3 hours, one needs at least 50 G. 
In comparison, leptonic modcls would givc a minimum jet luminosity of ^ 10^^ erg/s for 
this flare (Ghisellini 1998), about two orders of magnitude lower. This is duc to the much 
lower magnetic fields invoked there, and the less massive particles which drive the kinetic 
flow. 

In the framework of the jet-disk symbiosis (e.g. Falcke & Biermann 1995), the jet 
luminosity should not exceed the total accrction power Qaccr for the equilibrium state. 
Accretion theory relates the disk luminosity to the accretion power. Page & Thornc 
(1974) give L^i^k ~ (0.05 — 0.3)(5accr- Disk luminosities for 'typical' radio-loud AGN lie 
in the range L^isk ~ 10^ — 10^^ erg/s with BL Lac objects tending to the lower end 
on average. Specifically, for Mkn 501 there are no emission line measurements available, 
and this complicates the evaluation of the disk luminosity for this object. However, any 
observed UV-cmission in the flaring stage may put an upper limit on it. Again only 
historical data are available here, leaving room for speculation. We estimate L^isk ~ 
10^^ — 10^^ erg/s (Mufson et al 1984, Pian et al 1998), and obtain for the accretion power, 
Qaccr ~ (3 — 200) • 10^^ erg/s, at least a factor 5 below the necessary value to comply with 
the constraint of the disk-jet symbiosis. Note, however, that the estimate of (5accr is based 
on archival non- flaring data from Mkn 501, and we could speculate that either the disk 
has pushed more energy into the jet during TeV-flaring, or that the flaring stage can not 
be considered as a steady state. Also, accretion theory might predict larger conversion 
efficiencies of the accretion power into disk radiation than actually might occur in BL Lac 
objects. 

It is also instructive to consider the radiative efficiency of the proposed model in 
comparison to alternative models. Using Eq. 29 wc estimate the radiation efficiency of 
the protons we find (p ~ 0.01. This is similar to the value of (p in the Proton-Blazar model 
proposed by Mannheim (1993), and is fuUy in agreement with the results presented by 
Celotti & Fabian (1993) on a basis of a sample of 105 sources. In comparison, leptonic 
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models would give two orders of magnitude higher radiation efRciencies. 
4.1 Neutrino spectra 

Unlike leptonic models, hadronic blazar models may result in neutrino emission through 
the production and decay of charged mesons, e.g. tt'^ ^ /i'^ + foUowed by — > 

+ fe + T^^- Neutrinos escape without further interaction, and the predicted neutrino 
spectrum from Mkn 501 during flaring is shown in Fig. 8. We calculate the z/-emission from 
the source itself, and do not include here any additional contribution from escaping cosmic 
rays interacting while propagating through the cosmic microwave background radiation 
(see e.g. Protheroe and Johnson 1996). 

The proton injection spectrum is modified by the photo-hadronic interaction ratc 
which approximately follows a 7^'^'^ power-law for proton energies above ~ 10^ GeV (scc 
Fig. 3), where the nucleons interact preferably in the fiatter part of the target photon 
spectrum (e < 1.6 ke V) to produce mesons. This causes the resulting neutrino spectrum 
above Ej, 10^ GeV to be dN/dE,^ oc -B"^'^, whereas below this energy, the target 
photon field for pion production is the steeper part (1.6 keV < e < 42 ke V), and the 
corresponding neutrino spectrum is dN/dEy oc E~^''^. At even lower energies a further 
fiattening is due to pion production by the lowest energy protons at threshold. There is 
a steepening in the spectrum above 10^ GeV which needs some comment here. Neutrinos 
with energy < 10^ GeV are mostly produced near threshold in the A-resonance region, 
while the higher energy neutrinos are mainly produced in the secondary resonant region 
of the pion production cross section. 7r~ production is suppressed in the A-resonance, and 
hence we find a suppression of i/g-emission in comparison to i^e emission at low energies, 
whereas tt"*" and 7r~ multiplicities are comparable at higher CM energies, leading to roughly 
equal fluxes of and Pe- The effects of //^-synchrotron emission show up as a break at 
~ 10^ GcV in the observer's frame, whereas 7r^-synchrotron emission turns out to be 
unimportant in our model. 

Another important source of high energy neutrinos is the production and decay of 
charged kaons. In the case of p — 7-interactions positively charged kaons are produced, 
and are important in the secondary resonance region of the cross section (Miicke et al 
1999, Miicke et al 2000). They decay in ~ 64% of all cases into muons and dircct high 
energy muon-neutrinos. In contrast to the neutrinos originating from vr^ and //"^-decay, 
these muon-neutrinos will not have suffered energy losses through tt^- and /^^-synchrotron 
radiation, and therefore appear as an excess in comparison to the remaining neutrino 
flavors at the high energy end (E^, f« 10^-10^° GeV) of the emerging neutrino spectrum, 
and in addition cause the total neutrino spectrum to extend to ~ 10^° GeV. 

In contrast to previous SPB-jet models in which one expected equal photon and neu- 
trino energy fiuxes (e.g. Mannheim 1993, 1995), our model predicts a peak neutrino 
energy flux approximately two orders of magnitude lower than high energy gamma-rays. 
This is due to the synchrotron losses of protons being dominant, leading to gamma-ray 
emission at the expense of neutrino emission. 



17 



5 Conclusions 



This paper describes an application of our newly-developed Monte Carlo program which 
simulates a modified version of the stationary SPB-model. The Monte Carlo technique 
allows us to use exact cross sections, and all important emission processes are considered 
hcrc. As an example, we have used our codc to model the giant April 1997 TeV-flare 
from Mkn 501. Here, the TeV-photons are due to synchrotron radiation of the relativistic 
protons in the highly magnetized emission region. This proton synchrotron model was 
first proposed by Miicke & Protheroe (1999); a similar model for TeV emission in Mrk 501 
has just been proposed by Aharonian (2000), and his conclusions regarding the required 
Dopplcr factor and magnetic ficld are very similar to ours. Our model departs from 
the Standard SPB model as introduced first by Mannheim and co-workers mainly in two 
areas: (i) we use the observed synchrotron radiation as the target photon field for py- 
interactions and pair-synchrotron cascades assuming it to be produced by co-accelerated 
electrons, and (ii) our model takes into account synchrotron radiation from muons and 
protons. The model parameters derived assuming diffusive shock acceleration of e~ and 
p in a Kolmogorov turbulence spectrum are consistent with the X-ray to TeV-data in 
the flare state. However, the total jet power we obtain is too large to comply with the 
steady-state jet-disk symbiosis scenario, but then we are not dealing with a stready-state 
phenomenon. 

While the emerging cascade SED initiated by 7r° decay and n"^ synchrotron photons 
turns out to be relatively featureless, as was also found by, e. g., Mannheim (1993), the 
(see also Rachen 1999, and Rachen & Mannheim 2000) and, more importantly, the proton 
synchrotron radiation and its cascade produces a double-humped SED as is commonly 
observed in flaring blazars. For the present model, we find that proton synchrotron 
radiation dominates the TeV emission, while the contribution of the synchrotron radiation 
from the pairs, produced by photon-photon interactions of gamma-rays from the high 
energy hump is only minor. Our model considers the emission region to be homogcneous. 
Inhomogeneities in particle density, magnetic field, etc. within the source would result in a 
broader X-ray and TeV-peak in the SED. This indicates that for Mkn 501 a homogeneous 
model of the emission region seems to be appropriate. 

Being a hadronic model, our model predicts neutrino emission and we give the ex- 
pected neutrino fiux of Mrk 501 during flaring. Comparing our predicted neutrino flux 
with that for previous proton blazar models (e. g. Mannheim 1993), we find that the neu- 
trino output in our model is significantly less than was previously estimated due to the 
synchrotron losses dominating the energy losses of protons, producing synchrotron 7-rays 
at the expense of 7r° 7-rays and neutrinos. 
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Figure 1: (a) y F(9i,7jp^raax) versus shock angle (dashed lines) for different r]p values 
(numbers attached to the curves) and Tc = 4. Note that the diffusion approximation limits 
the mean free path to r) < Region below the solid line satisfies iyar-D > iocc,p,max 

for a typical set of TeV-blazar parameters B ^ 20 G, ~ 10, ui = 0.5c, = 1, rc = 4, 
tyar = 12 hours. (b) Ratio F{6i,rje,^ax) / F{6i, 7yp,max) versus 6i for different turbulence 
spectra. 
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Figure 2: AUowed parameter space (shaded area) for 7p,max, 7e,max in the SPB-model for 
different magnetic turbulence spetra. The diamond symbol corresponds to the Mkn 501- 
model presented below {B Ri 20 G, D ^ 10, Ui — 0.5c, P — 1, ^var — 12 hours). 
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Figure 3: Mean energy loss time of p for 7r-photoproduction (tt), Bethe-Heitler pair 
production (e"'"e~) and synchrotron radiation (syn). Loss times for tt^- and /i"^ for syn- 
chrotron radiation (syn tt, syn /i) are also shown and compared with their mean decay 
time scales (decay tt, decay /i). The acceleration time scale (acc), based on Kolmogorov 
turbulence, is calculated for = 4, Ui = 0.5c and shock angle 6i = 85°. The curvature 
of the acceleration time scale is caused by the non-linear dependence of F[6, rjp) on 7. 
The adiabatic loss time (adiab) is assumed to be 2\B/IS\ ^ R/ui -Dtvar- We adopt 
B 20 G, and all quantities are in the co-moving frame of the jet. 
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Figure 4: (a) Average cascade spectrum initiated by 7r°-decay, and (b) tt"*^ — > /^"^ — >• e^- 
decay synchrotron photons (lower panel). Dot-dashed and dotted histogram shows the 
Ist and 2nd generation synchrotron photons on production. Sohd histogram shows the 
emerging photon spectrum on production. 
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Figure 5: (a) Average emerging cascade spectrum initiated by Bethe-Heitler pair pro- 
duction. (b) Average emerging cascade spectrum initiated by p- (solid histogram) and 
//='=-synchrotron photons (dot-dashed histogram) . Note the different vertical scale in both 
figures. 
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Figure 6: Best-fit model (histogram) in comparison with the data of the 16 April 1997- 
fiare of Mkn 501. Photon absorption on the cosmic diffuse background radiation field is 

not includcd. Straight solid lincs: parametrization of the obscrved synchrotron spcctrum 
(BcppoSAX & OSSE) and obscrved TeV-emission corrected for cosmic background ab- 
sorption (Bednarek & Protheroe 1999); the 100 MeV upper limit is from Catanese et al 
1997 (observed 9-15 April 1997); dashed-dotted hne: input target spectrum. 
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Figurc 7: Dependence of total jet luminosity Ljet = LB+Lp^hot+-^>p,coid+-^c and its different 
contributions, L b for the magnetic field, Lp^hot and Lp^coid for the hot and cold protons 
and Le for the electron contribution, on the Doppler factor D for model parameters which 
reasonably fit the Mkn 501 flare SED. Fixed parameters are tyar = 12 hours, /3i = 0.5 
and Np oc 7"^ for 2 < 7^ < 7p,max = 3 x 10^*^, adjusted parameters are B, R, Up, «target) 
the synchrotron photon break energy 65 and its maximum photon energy 62- Ljet shows 
a minimum of ~ lO^^erg/s at around D w 12. L40 — L/ 10^° erg/s. 
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Figure 8: Predicted spectra of z/g, ^'e, and from Mkn 501. The contribution of vs 
due to pion production of the emerging cosmic rays whilc propagating through the cosmic 
background is not considered here. 
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